/*
 *  glfem.c
 *  Library for MECA1120 : Finite Elements for dummies
 *
 *  Copyright (C) 2012 UCL-IMMC : Vincent Legat
 *  All rights reserved.
 *
 */

#include "glfem.h"

glfemPlot thePlot;

void GLFWCALL glfemPlotReshape(int w, int h);
void glMakeRasterFont(void);
void setColor(float value, int numberOfColors);
void getColor(double value, int numberOfColors, float* R, float* G, float* B);
double scale(double minimum, double maximum, double value);


static int gRasterH = 800;
static int gRasterV = 600;
static int numberColors = 50;

GLubyte space[] = 
    {0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00};

GLubyte letters[][13] = {
    {0x00, 0x00, 0xc3, 0xc3, 0xc3, 0xc3, 0xff, 0xc3, 0xc3, 0xc3, 0x66, 0x3c, 0x18}, // A
    {0x00, 0x00, 0xfe, 0xc7, 0xc3, 0xc3, 0xc7, 0xfe, 0xc7, 0xc3, 0xc3, 0xc7, 0xfe}, // B
    {0x00, 0x00, 0x7e, 0xe7, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xe7, 0x7e}, // C
    {0x00, 0x00, 0xfc, 0xce, 0xc7, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc7, 0xce, 0xfc}, // D
    {0x00, 0x00, 0xff, 0xc0, 0xc0, 0xc0, 0xc0, 0xfc, 0xc0, 0xc0, 0xc0, 0xc0, 0xff}, // E
    {0x00, 0x00, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xfc, 0xc0, 0xc0, 0xc0, 0xff}, // F
    {0x00, 0x00, 0x7e, 0xe7, 0xc3, 0xc3, 0xcf, 0xc0, 0xc0, 0xc0, 0xc0, 0xe7, 0x7e}, // G
    {0x00, 0x00, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xff, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3}, // H
    {0x00, 0x00, 0x7e, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x7e}, // I
    {0x00, 0x00, 0x7c, 0xee, 0xc6, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06}, // J
    {0x00, 0x00, 0xc3, 0xc6, 0xcc, 0xd8, 0xf0, 0xe0, 0xf0, 0xd8, 0xcc, 0xc6, 0xc3}, // K
    {0x00, 0x00, 0xff, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0}, // L
    {0x00, 0x00, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xdb, 0xff, 0xff, 0xe7, 0xc3}, // M
    {0x00, 0x00, 0xc7, 0xc7, 0xcf, 0xcf, 0xdf, 0xdb, 0xfb, 0xf3, 0xf3, 0xe3, 0xe3}, // N
    {0x00, 0x00, 0x7e, 0xe7, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xe7, 0x7e}, // O
    {0x00, 0x00, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0xfe, 0xc7, 0xc3, 0xc3, 0xc7, 0xfe}, // P
    {0x00, 0x00, 0x3f, 0x6e, 0xdf, 0xdb, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x66, 0x3c}, // Q
    {0x00, 0x00, 0xc3, 0xc6, 0xcc, 0xd8, 0xf0, 0xfe, 0xc7, 0xc3, 0xc3, 0xc7, 0xfe}, // R
    {0x00, 0x00, 0x7e, 0xe7, 0x03, 0x03, 0x07, 0x7e, 0xe0, 0xc0, 0xc0, 0xe7, 0x7e}, // S
    {0x00, 0x00, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0xff}, // T
    {0x00, 0x00, 0x7e, 0xe7, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3}, // U
    {0x00, 0x00, 0x18, 0x3c, 0x3c, 0x66, 0x66, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3}, // V
    {0x00, 0x00, 0xc3, 0xe7, 0xff, 0xff, 0xdb, 0xdb, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3}, // W
    {0x00, 0x00, 0xc3, 0x66, 0x66, 0x3c, 0x3c, 0x18, 0x3c, 0x3c, 0x66, 0x66, 0xc3}, // X
    {0x00, 0x00, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x3c, 0x3c, 0x66, 0x66, 0xc3}, // Y
    {0x00, 0x00, 0xff, 0xc0, 0xc0, 0x60, 0x30, 0x7e, 0x0c, 0x06, 0x03, 0x03, 0xff}  // Z
};

GLubyte lowletters[][13] = {
    {0x00, 0x00, 0x7d, 0xc3, 0xc3, 0xc3, 0x7f, 0x03, 0x7e, 0x00, 0x00, 0x00, 0x00}, // a
    {0x00, 0x00, 0xfe, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xfe, 0xc0, 0xc0, 0xc0, 0xc0}, // b
    {0x00, 0x00, 0x7f, 0xc0, 0xc0, 0xc0, 0xc0, 0xc0, 0x7f, 0x00, 0x00, 0x00, 0x00}, // c
    {0x00, 0x00, 0x7f, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x7f, 0x03, 0x03, 0x03, 0x03}, // d
    {0x00, 0x00, 0x7e, 0xc0, 0xc0, 0xfe, 0xc3, 0xc3, 0x7e, 0x00, 0x00, 0x00, 0x00}, // e
    {0x00, 0x00, 0x3c, 0x18, 0x18, 0x18, 0x18, 0x18, 0x7e, 0x18, 0x18, 0x18, 0x0e}, // f
    {0x7f, 0x03, 0x7f, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x7e, 0x00, 0x00, 0x00, 0x00}, // g
    {0x00, 0x00, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xfe, 0xc0, 0xc0, 0xc0, 0xc0}, // h
    {0x00, 0x00, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x00, 0x18, 0x18, 0x00}, // i
    {0x70, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x00, 0x18, 0x18, 0x00}, // j
    {0x00, 0x00, 0xc3, 0xc7, 0xce, 0xfc, 0xfe, 0xc7, 0xc3, 0xc0, 0xc0, 0xc0, 0xc0}, // k
    {0x00, 0x00, 0x0c, 0x1c, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18}, // l
    {0x00, 0x00, 0xdb, 0xdb, 0xdb, 0xdb, 0xdb, 0xdb, 0xfe, 0x00, 0x00, 0x00, 0x00}, // m
    {0x00, 0x00, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xfe, 0x00, 0x00, 0x00, 0x00}, // n
    {0x00, 0x00, 0x7e, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x7e, 0x00, 0x00, 0x00, 0x00}, // o
    {0xc0, 0xc0, 0xfe, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xfe, 0x00, 0x00, 0x00, 0x00}, // p
    {0x03, 0x03, 0x7f, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x7f, 0x00, 0x00, 0x00, 0x00}, // q
    {0x00, 0x00, 0xc0, 0xc0, 0xc0, 0xc0, 0xe0, 0xf0, 0xdf, 0x00, 0x00, 0x00, 0x00}, // r
    {0x00, 0x00, 0xfe, 0x03, 0x03, 0x7e, 0xc0, 0xc0, 0x7f, 0x00, 0x00, 0x00, 0x00}, // s
    {0x00, 0x00, 0x0e, 0x18, 0x18, 0x18, 0x18, 0x18, 0x7e, 0x18, 0x18, 0x18, 0x18}, // t
    {0x00, 0x00, 0x7f, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x00, 0x00, 0x00, 0x00}, // u
    {0x00, 0x00, 0x18, 0x3c, 0x66, 0x66, 0xc3, 0xc3, 0xc3, 0x00, 0x00, 0x00, 0x00}, // v
    {0x00, 0x00, 0x66, 0x7e, 0xdb, 0xdb, 0xdb, 0xdb, 0xdb, 0x00, 0x00, 0x00, 0x00}, // w
    {0x00, 0x00, 0xc3, 0xe7, 0x3c, 0x18, 0x3c, 0xe7, 0xc3, 0x00, 0x00, 0x00, 0x00}, // x
    {0x7f, 0x03, 0x7f, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0xc3, 0x00, 0x00, 0x00, 0x00}, // y
    {0x00, 0x00, 0xff, 0xc0, 0x70, 0x1c, 0x06, 0x03, 0xff, 0x00, 0x00, 0x00, 0x00}, // z
};

GLubyte numletters[][13] = {
    {0x00, 0x00, 0x3c, 0x66, 0x66, 0x66, 0x66, 0x66, 0x66, 0x66, 0x66, 0x66, 0x3c}, // 0
    {0x00, 0x00, 0x3c, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x78, 0x38, 0x18}, // 1
    {0x00, 0x00, 0x7e, 0x60, 0x60, 0x60, 0x60, 0x3c, 0x06, 0x06, 0x66, 0x66, 0x3c}, // 2
    {0x00, 0x00, 0x3c, 0x66, 0x06, 0x06, 0x06, 0x1c, 0x06, 0x06, 0x06, 0x66, 0x3c}, // 3
    {0x00, 0x00, 0x06, 0x06, 0x06, 0x06, 0x06, 0x7f, 0x66, 0x36, 0x1e, 0x0e, 0x06}, // 4
    {0x00, 0x00, 0x3c, 0x66, 0x06, 0x06, 0x06, 0x7c, 0x60, 0x60, 0x60, 0x60, 0x7e}, // 5
    {0x00, 0x00, 0x3c, 0x66, 0x66, 0x66, 0x66, 0x66, 0x7c, 0x60, 0x60, 0x66, 0x3c}, // 6
    {0x00, 0x00, 0x0c, 0x0c, 0x0c, 0x0c, 0x0c, 0x1f, 0x06, 0x06, 0x06, 0x06, 0x7e}, // 7
    {0x00, 0x00, 0x3c, 0x66, 0x66, 0x66, 0x66, 0x3c, 0x66, 0x66, 0x66, 0x66, 0x3c}, // 8
    {0x00, 0x00, 0x3c, 0x66, 0x06, 0x06, 0x06, 0x3e, 0x66, 0x66, 0x66, 0x66, 0x3c}, // 9
    {0x00, 0x00, 0x00, 0x18, 0x18, 0x00, 0x00, 0x18, 0x18, 0x00, 0x00, 0x00, 0x00}, // :
    {0x00, 0x00, 0x30, 0x18, 0x18, 0x00, 0x00, 0x18, 0x18, 0x00, 0x00, 0x00, 0x00}, // ;
    {0x00, 0x00, 0x06, 0x1c, 0x30, 0x60, 0x30, 0x1c, 0x06, 0x00, 0x00, 0x00, 0x00}, // <
    {0x00, 0x00, 0x00, 0x00, 0x3c, 0x00, 0x00, 0x3c, 0x00, 0x00, 0x00, 0x00, 0x00}, // =
    {0x00, 0x00, 0x60, 0x38, 0x0c, 0x06, 0x0c, 0x38, 0x60, 0x00, 0x00, 0x00, 0x00}, // >
    {0x00, 0x00, 0x18, 0x18, 0x00, 0x18, 0x18, 0x18, 0x0c, 0x06, 0x06, 0x66, 0x3c}, // ?
};

GLubyte specialletters[][13] = {
    {0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00}, // space
    {0x00, 0x00, 0x18, 0x18, 0x00, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18, 0x18}, // !
    {0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x24, 0x24, 0x00, 0x00}, // "
    {0x00, 0x00, 0x24, 0x24, 0x7e, 0x7e, 0x24, 0x7e, 0x7e, 0x24, 0x24, 0x00, 0x00}, // #
    {0x00, 0x00, 0x18, 0x3c, 0x5a, 0x5a, 0x1a, 0x3c, 0x58, 0x58, 0x5a, 0x3c, 0x18}, // $
    {0x00, 0x00, 0x44, 0x4a, 0x6a, 0x24, 0x30, 0x18, 0x0c, 0x24, 0x56, 0x52, 0x22}, // %
    
    {0x00, 0x00, 0x79, 0xcf, 0xc6, 0xcf, 0x79, 0x70, 0x78, 0xcc, 0xcc, 0xcc, 0x78}, // &
    {0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x08, 0x08, 0x18, 0x00, 0x00}, // '
    {0x00, 0x00, 0x0c, 0x18, 0x18, 0x30, 0x30, 0x30, 0x30, 0x30, 0x18, 0x18, 0x0c}, // (
    {0x00, 0x00, 0x30, 0x18, 0x18, 0x0c, 0x0c, 0x0c, 0x0c, 0x0c, 0x18, 0x18, 0x30}, // )
    {0x00, 0x00, 0x00, 0x00, 0x10, 0x54, 0x38, 0x54, 0x10, 0x00, 0x00, 0x00, 0x00}, // *
    {0x00, 0x00, 0x00, 0x00, 0x10, 0x10, 0x7c, 0x10, 0x10, 0x00, 0x00, 0x00, 0x00}, // +
    {0x00, 0x30, 0x18, 0x18, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00}, // ,
    {0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x3c, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00}, // -
    {0x00, 0x00, 0x18, 0x18, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00}, // .
    {0x00, 0x00, 0x60, 0x60, 0x30, 0x30, 0x18, 0x18, 0x18, 0x0c, 0x0c, 0x06, 0x06}, // /
};

GLuint fontOffset;

void glMakeRasterFont(void)
{
    GLuint i, j;
    glPixelStorei(GL_UNPACK_ALIGNMENT, 1);  
    glShadeModel (GL_FLAT);
    fontOffset = glGenLists (128);
    
    for (i = 0,j = 'A'; i < 26; i++,j++) {
        glNewList(fontOffset + j, GL_COMPILE);
        glBitmap(8, 13, 0.0, 2.0, 10.0, 0.0, letters[i]);
        glEndList(); }
    
    for (i = 0,j = 'a'; i < 26; i++,j++) {
        glNewList(fontOffset + j, GL_COMPILE);
        glBitmap(8, 13, 0.0, 2.0, 10.0, 0.0, lowletters[i]);
        glEndList(); }
    
    for (i = 0,j = '0'; i < 16; i++,j++) {
        glNewList(fontOffset + j, GL_COMPILE);
        glBitmap(8, 13, 0.0, 2.0, 10.0, 0.0, numletters[i]);
        glEndList(); }
    
    for (i = 0,j = ' '; i < 16; i++,j++) {
        glNewList(fontOffset + j, GL_COMPILE);
        glBitmap(8, 13, 0.0, 2.0, 10.0, 0.0, specialletters[i]);
        glEndList(); }
    glShadeModel (GL_SMOOTH);
}

void glfemSetRasterSize(int h, int v)
{
    gRasterH = h;
    gRasterV = v;
}

void glfemDrawMessage(int h, int v, char *s)
{
    // Les coordonnées négatives sont normalement admises :-)
    // On conserve le strint entier même si le début est off-screen

    int off;    
    glPushAttrib(GL_ALL_ATTRIB_BITS);   
    glDisable(GL_TEXTURE_2D);   
    glMatrixMode(GL_PROJECTION);
    glPushMatrix();
    glLoadIdentity();
    glOrtho (0.0, gRasterH, gRasterV, 0.0, -1.0, 1.0);
    glMatrixMode(GL_MODELVIEW);
    glPushMatrix();
    glLoadIdentity();
    
    if (h < 0)   h = gRasterH + h - strlen(s)*10;
    if (v < 0)   v = gRasterV + v; 
    glRasterPos2i(h, v);    
    glListBase(fontOffset);
    
    if (h >= 0) {
      glRasterPos2i(h, v);
      glCallLists((GLsizei) strlen(s), GL_UNSIGNED_BYTE, (GLubyte *) s); }
    else { 
      off = (h-9)/10;
      glRasterPos2i(h - off*10, v);
      if (strlen(s)+off > 0) glCallLists((GLsizei) strlen(s)+off, GL_UNSIGNED_BYTE, (GLubyte *) &s[-off]); }
    
    glPopMatrix(); 
    glMatrixMode(GL_PROJECTION);
    glPopMatrix();
    glMatrixMode(GL_MODELVIEW);
    glPopAttrib ();
}

void setColor(float value, int numberOfColors)
{
    float R,G,B;
    getColor(value, numberOfColors, &R, &G, &B);
    glColor3f(R,G,B);
}

void getColor(double value, int numberOfColors, float* R, float* G, float* B)
{
    if (value > 1) value = 1;
    value = value * (numberOfColors);
    value = (int)(value - 0.00000001);
    value = value / (numberOfColors - 1);
    
    value = 1 - value;
    if (value < 0) value=0;
    if (value > 1) value=1; 
    *R = 3.5*(1-value)*(1-value);
    *G = (1-value)*(value)*3.5;
    *B = value*value;
}

double scale(double minimum, double maximum, double value)
{
    if (value < minimum)        return 0;
    if (minimum == maximum)     return minimum;
    return (value - minimum) / fabs(maximum - minimum);
}

void glfemDrawNodes(double* x, double* y,int n) 
{
    int i;
    glEnable(GL_POINT_SMOOTH);
    glPointSize(12.0);
    glBegin(GL_POINTS);
    for (i = 0; i < n; i++) {      
        glVertex2f(x[i],y[i]); }
    glEnd();
}

void glfemDrawCurve(femVector *theX, femVector *theY)
{
    int i;
    int size = theX->size;    
    glBegin(GL_LINES);
    for (i = 0; i < size-1; i++) {
        glVertex2d(theX->X[i], theY->X[i]); 
        glVertex2d(theX->X[i+1], theY->X[i+1]);}
    glEnd();

}

void glfemDrawCurveDiscrete(femVector *theX, femVector *theY)
{
    glfemDrawCurve(theX, theY);
    glfemDrawNodes(theX->X,theY->X,theX->size);
    
}


double interpolate(double u[3], double xsi, double eta)
{
	return u[0]*xsi + u[1]*eta + u[2]*(1.0-xsi-eta);
}

void glfemDrawColorRecursiveTriangle(double *x, double *y, 
									 double *xsi, double *eta, 
									 double (*f)(int,double,double,double,double),
									 double fMin, double fMax,
									 int iElem,
									 int n)
{
    int i,j;
    const int nodes[4][3] = {{0,3,5},{3,1,4},{5,4,2},{3,4,5}};
    const double xsiRef[6] = {0.0,1.0,0.0,0.5,0.5,0.0};
    const double etaRef[6] = {0.0,0.0,1.0,0.0,0.5,0.5};
    double xSub[3], xsiSub[3];
    double ySub[3], etaSub[3];
	float coord[6];
	double values[3];
    
    if (n <= 0) {
        for (j=0; j<3; j++) {
			values[j] = scale(fMin,fMax,f(iElem,xsi[j],eta[j],x[j],y[j])); 
			coord[j*2+0] = x[j];
            coord[j*2+1] = y[j]; }                 
		glfemDrawColorTriangle(coord,values);  }
    else {
		for (i=0; i<4; i++) {
			for (j=0; j<3; j++) {
				double xsiLoc = xsiRef[nodes[i][j]];
				double etaLoc = etaRef[nodes[i][j]];
				xSub[j] = interpolate(x,xsiLoc,etaLoc);
				ySub[j] = interpolate(y,xsiLoc,etaLoc);  
				xsiSub[j] = interpolate(xsi,xsiLoc,etaLoc);
				etaSub[j] = interpolate(eta,xsiLoc,etaLoc); } 
			glfemDrawColorRecursiveTriangle(xSub,ySub,xsiSub,etaSub,f,fMin,fMax,iElem,n-1); }}
}

void glfemDrawColorTriangle(float* coord, double *values) 
{   
    GLfloat colors[9];
    getColor(values[0],numberColors,&colors[0],&colors[1],&colors[2]);
    getColor(values[1],numberColors,&colors[3],&colors[4],&colors[5]);
    getColor(values[2],numberColors,&colors[6],&colors[7],&colors[8]);
    glEnableClientState(GL_VERTEX_ARRAY);
    glEnableClientState(GL_COLOR_ARRAY);
    glVertexPointer(2, GL_FLOAT, 0, coord);
    glColorPointer(3, GL_FLOAT, 0, colors);
    glDrawArrays(GL_TRIANGLES, 0, 3);
    glDisableClientState(GL_COLOR_ARRAY);
    glDisableClientState(GL_VERTEX_ARRAY);
}

void glfemDrawTriangle(float* coord) 
{
    glEnableClientState(GL_VERTEX_ARRAY);
    glVertexPointer(2, GL_FLOAT, 0, coord);
    glDrawArrays(GL_LINE_LOOP, 0, 3);
    glDisableClientState(GL_VERTEX_ARRAY);
}

void glfemReshapeWindows(int w, int h)
{
    double minX  = thePlot.Xmin;
    double maxX  = thePlot.Xmax;
    double minY  = thePlot.Ymin;
    double maxY  = thePlot.Ymax;
    double sizeX = (maxX-minX)/1.9;
    double meanX = (maxX+minX)/2.0; 
    double sizeY = (maxY-minY)/1.9;
    double meanY = (maxY+minY)/2.0;
    
    double ratio = (GLfloat) h / (GLfloat) w;
    double size = Max(sizeX,sizeY);
    double left = meanX - size;
    double right = meanX + size;
    double bottom = meanY- size*ratio;
    double top = meanY + size*ratio;        
        
    glViewport(0,0,(GLsizei) w, (GLsizei) h);    
    glClearColor( 0.9f, 0.9f, 0.8f, 0.0f );
    //glClearColor( 1.0f, 1.0f, 1.0f, 0.0f );
    glClear(GL_COLOR_BUFFER_BIT);   
    glMatrixMode(GL_PROJECTION);
    glLoadIdentity();
    glOrtho((GLdouble)left, (GLdouble)right, (GLdouble)bottom, (GLdouble)top, -5.0, 5.0);
    glMatrixMode(GL_MODELVIEW);
    glLoadIdentity();
}


void glfemPlotField(femMesh *theMesh, double *u)
{
    int i,j,*nodes;
    float coord[6];
    double uLoc[3];
    double uMax = femMax(u,theMesh->nNode);
    double uMin = femMin(u,theMesh->nNode);
    
    for (i = 0; i < theMesh->nElem; ++i) {
        nodes = femMeshElem(theMesh,i);
        for (j=0; j < 3; ++j) {
            coord[j*2+0] = theMesh->X[nodes[j]-1];
            coord[j*2+1] = theMesh->Y[nodes[j]-1];
            uLoc[j] = scale(uMin,uMax,u[nodes[j]-1]); 
		}              
        glfemDrawColorTriangle(coord,uLoc); }
    glColor3f(0.0, 0.0, 0.0);
    glfemPlotMesh(theMesh);
}

void glfemPlotDiscField(femMesh *theMesh, double (*f)(int,double,double,double,double),double fMin, double fMax,int level)
{
    int i,j,*nodes;
    double x[3],y[3];
	double xsi[3] = {0.0,1.0,0.0};
    double eta[3] = {0.0,0.0,1.0};
    
    for (i = 0; i < theMesh->nElem; ++i) {
        nodes = femMeshElem(theMesh,i);
        for (j=0; j < 3; ++j) {
            x[j] = theMesh->X[nodes[j]-1];
            y[j] = theMesh->Y[nodes[j]-1]; }              
        glfemDrawColorRecursiveTriangle(x,y,xsi,eta,f,fMin,fMax,i,level); }
    glColor3f(0.0, 0.0, 0.0);
    glfemPlotMesh(theMesh);
}


void glfemPlotMesh(femMesh *theMesh)
{
    int i,j,*nodes;
    float coord[6];
    
    for (i = 0; i < theMesh->nElem; ++i) {
        nodes = femMeshElem(theMesh,i);
        for (j=0; j < 3; ++j) {
            coord[j*2+0] = theMesh->X[nodes[j]-1];
            coord[j*2+1] = theMesh->Y[nodes[j]-1]; }                 
        glfemDrawTriangle(coord); }
}

void glfemPlotEdges(femEdges *theEdges)
{
    int i;
    femMesh* theMesh = theEdges->mesh;
    glBegin(GL_LINES);
    for (i = 0; i < theEdges->nEdge; i++) {
        glVertex2d(theMesh->X[theEdges->edges[i].node[0] - 1], theMesh->Y[theEdges->edges[i].node[0] - 1]);
        glVertex2d(theMesh->X[theEdges->edges[i].node[1] - 1], theMesh->Y[theEdges->edges[i].node[1] - 1]); }
    glEnd();
}

void glfemPlotBnd(femEdges *theEdges)
{
    int i;
    femMesh* theMesh = theEdges->mesh;
    glBegin(GL_LINES);
    for (i = 0; i < theEdges->nEdge; i++) {
        if (theEdges->edges[i].elem[1] < 0) {
            glVertex2d(theMesh->X[theEdges->edges[i].node[0] - 1], theMesh->Y[theEdges->edges[i].node[0] - 1]);
            glVertex2d(theMesh->X[theEdges->edges[i].node[1] - 1], theMesh->Y[theEdges->edges[i].node[1] - 1]); }}
    glEnd();
}


void glfemMessage(char *message)
{
  //  glColor3f (0.0, 0.0, 0.0);
    glfemDrawMessage(20,40,message);   
 //   glfemDrawMessage(20,460,message);   
}

void glfemInit(char *theWindowName, femMesh *theMesh, femEdges *theEdges, double *theField, char *theMessage)
{
    thePlot.mesh = theMesh;
    thePlot.edges = theEdges;
    thePlot.field = theField;
    thePlot.message = theMessage;
    glfwInit();
    glfwOpenWindow(480,480,0,0,0,0,0,0,GLFW_WINDOW);    
    glfemSetRasterSize(480,480);
    glfwSetWindowTitle(theWindowName);
    glShadeModel(GL_SMOOTH);
    glMakeRasterFont();
    glfemSetViewPortMesh(theMesh);
}

void glfemSetDrawFunction(void (*draw)(void))
{
    thePlot.draw = draw;
}

void glfemSetAnimateFunction(void (*animate)(void))
{
    thePlot.animate = animate;
}

void glfemSetResetFunction(void (*reset)(void))
{
    thePlot.reset = reset;
}

void glfemSetKeyFunction(void (*key)(char))
{
    thePlot.key = key;
}


void glfemSetViewPort(double Xmin, double Xmax, double Ymin, double Ymax)
{
    thePlot.Xmin = Xmin;
    thePlot.Xmax = Xmax;
    thePlot.Ymin = Ymin;
    thePlot.Ymax = Ymax;
}

void glfemSetViewPortMesh(femMesh* theMesh)
{
    thePlot.mesh = theMesh;
    thePlot.Xmin = femMin(theMesh->X,theMesh->nNode);
    thePlot.Xmax = femMax(theMesh->X,theMesh->nNode);
    thePlot.Ymin = femMin(theMesh->Y,theMesh->nNode);
    thePlot.Ymax = femMax(theMesh->Y,theMesh->nNode);
}


void glfemMainLoop(char *theWindowName)
{	
    glfwInit();
    glfwOpenWindow(480,480,0,0,0,0,0,0,GLFW_WINDOW);    
    glfemSetRasterSize(480,480);
    glfwSetWindowTitle(theWindowName);
    glShadeModel(GL_SMOOTH);
    glMakeRasterFont();
     do { //glfwWaitEvents();
        
        int width, height;
        glfwGetWindowSize( &width, &height );
		
        glfemReshapeWindows(width,height);
    //   glfwSleep(10);
        thePlot.draw();         
        thePlot.animate();         
        glfwSwapBuffers();

        if (glfwGetKey('R') == GLFW_PRESS)    thePlot.key('r');
        if (glfwGetKey('U') == GLFW_PRESS)    thePlot.key('u');
        if (glfwGetKey('V') == GLFW_PRESS)    thePlot.key('v');
        if (glfwGetKey('E') == GLFW_PRESS)    thePlot.key('e');
        if (glfwGetKey('S') == GLFW_PRESS)    thePlot.key('s');
        if (glfwGetKey('N') == GLFW_PRESS)    thePlot.key('n');
		if (glfwGetKey(GLFW_KEY_ESC) == GLFW_PRESS) break;
		if (glfwGetWindowParam(GLFW_OPENED) == GL_FALSE) break;
    } while(1);
    
    glfwTerminate();    
    exit(EXIT_SUCCESS); 
}
